Density dependence of the symmetry energy and the nuclear equation of state: A 

Dynamical and Statistical model perspective 
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The density dependence of the symmetry energy in the equation of state of isospin asymmetric 
nuclear matter is of significant importance for studying the structure of systems as diverse as the 
neutron-rich nuclei and the neutron stars. A number of reactions using the dynamical and the 
statistical models of multifragmentation, and the experimental isoscaling observable, is studied to 
extract information on the density dependence of the symmetry energy. It is observed that the 
dynamical and the statistical model calculations give consistent results assuming the sequential 
decay effect in dynamical model to be small. A comparison with several other independent studies 
is also made to obtain important constraint on the form of the density dependence of the symmetry 
energy. The comparison rules out an extremely " stiff " and " soft " form of the density dependence 
of the symmetry energy with important implications for astrophysical and nuclear physics studies. 

PACS numbers: 21.30.Fe, 25.70.-z, 25.70.Lm, 25.70.Mn, 25.70.Pq 



I. INTRODUCTION 

The fundamental goal of nuclear physics is to under- 
stand the basic building blocks of nature - neutrons and 
protons - and the nature of interaction that binds them 
together into nuclear matter. Studying the nature of 
matter and the strength of nuclear interaction is key to 
understanding some of the fundamental problems such 
as, How are elements formed? How do stars explode into 



supernova? What kind of matter exists inside a neu- 
tron star? How are neutrons compressed inside a neutron 
star to density trillions of times greater than on earth ? 
What determines the density-pressure relation, the so- 
called equation of state? 

The key ingredient for constructing the nuclear equa- 
tion of state is the basic nucleon-nucleon interaction. Un- 
til now our understanding of the nucleon-nucleon inter- 
action has come from studying nuclear matter that is 
symmetric in isospin (neutron-to-proton ratio, N/Z rj 1) 
and matter found near normal nuclear density (p « 0.16 
/m~ 3 ). It is not known how far this understanding re- 
mains valid as one goes away from the normal nuclear 
density and symmetric nuclear matter. Various interac- 
tions used in " ab initio " microscopic calculations pre- 
dict different forms of the nuclear equation of state above 
and below the normal nuclear matter density .and away 
from the symmetric nuclear matter [l|, H, H, |j, [E @1 • As a 
result, the symmetry energy, which is the difference in en- 
ergy between the pure neutron matter and the symmetric 
nuclear matter, shows very different behavior above and 
below normal nuclear density Q (see Fig. 1). 

In general, two different forms of the density depen- 
dence of the symmetry energy have been predicted. One, 
where the symmetry energy increases monotonically with 
increasing density (" stiff " dependence) and the other, 
where the symmetry energy increases initially up to nor- 
mal nuclear density and then decreases at higher densities 
( " soft " dependence) . Constraining the form of the den- 
sity dependence of the symmetry energy is important not 




FIG. 1: (Color online) Symmetry energy as a function of den- 
sity predicted by microscopic " ab initio " calculations. The 
left panel shows the low-density region, while the right panel 
displays the high-density range. The figure is taken from Ref. 
1- 



only for a better understanding of the nucleon-nucleon 
interaction, and hence its extrapolation to the structure 
of neutron-rich nuclei 0, @, 0, [323], but also for deter- 
mining the structure of com pac t stellar objects such as 
neutron stars [H G3, G3, Q Uj, QJ, [T3| • For example, a 
" stiff " form of the density dependence of the symmetry 
energy is predicted to lead to a large neutron skin thick- 
ness compared to a " soft " dependence Qil [H, 111, H3| • 
Similarly, a " stiff " dependence of the symmetry energy 
can result in rapid cooling of a neutron star, and a larger 
neutron star radius, compared to a " soft " density depen- 
dence of the symmetry energy [2(J HJ, • The nuclear 
Equation Of State (EOS) is therefore a fundamental en- 
tity that determines the properties of systems as diverse 
as atomic nuclei and neutron stars, and the knowledge of 
which is of significant importance [lj| Qjl [23| . 

Experimentally, the best possible means of studying 
the nuclear equation of state at sub-normal nuclear den- 
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sity is through intermediate-energy heavy-ion reactions 
[2a . l27l ]. In this kind of reaction, an excited nucleus 
(the composite of the projectile and the target nucleus) 
expands to a sub-nuclear density and disintegrates into 
various light and heavy fragments in a process called 
multifragmentation. By studying the isotopic yield dis- 
tribution of these fragments one can extract important 
information about the symmetry energy and its density 
dependence. Current studies on the nuclear equation of 
state are limited to beams consisting of stable nuclei. It 
is hoped that in the future radioactive beam facilities 
such as,FAIR (GSI) [H, SPIRAL2 (GANIL) and FRIB 
(USA) [25J will provide tremendous opportunities for ex- 
ploring the nuclear EOS in regions never before studied 
(i.e., extreme isospin and away from normal nuclear den- 
sity). 

In this work, we have made an attempt to study 
the density dependence of the symmetry energy using 
two different theoretical approaches for studying multi- 
fragmentation, namely the dynamical and the statistical 
model approaches of multifragmentation. In section II, 
the isoscaling technique used to study the density de- 
pendence of the symmetry energy, and their different in- 
terpretations in terms of statistical and dynamical ap- 
proach, are presented. In section III and IV, a brief 
description of the experiment and the experimental re- 
sults are presented. The dynamical and the statistical 
approaches used to interpret the experimental results are 
presented in section V. A comparsion between the two 
approaches with other independent studies is presented 
in section VI. Finally, a discussion and summary, and 
conclusions are presented in section VII and VIII, re- 
spectively. 

II. SYMMETRY ENERGY AND THE ISOTOPIC 
YIELD DISTRIBUTION 

It has been shown from experimental measurements 
that the ratio of the fragment yields, K 2 i(N,Z), taken 
from two different multifragmentation reactions, 1 and 
2, obeys an exponential dependence on the neutron num- 
ber (N) and the proton number (Z) of the fragments; an 
observation known as isoscaling [23, HH [3(| ■ The depen- 
dence is characterized by the relation, 

R 2 i(N,Z) = Y 2 (N,Z)/Y 1 (N,Z) = C e {aN+f)z ^ (1) 

Where, Y 2 and Y\ are the fragment yields from the 
neutron-rich and the neutron-deficient systems, respec- 
tively. C is an overall normalization factor, and a and /3 
are the parameters characterizing the isoscaling behavior. 
Isoscaling is also theoretically predicted by the dynam- 
ical [HI Bl [H, [H, [H and statistical [!& [13, \M, HI 
models of multifragmentation. In these models, the dif- 
ference in the chemical potential of systems with different 
neutron-to-proton ration {N/Z) is directly related to the 
isoscaling parameter a. The isoscaling parameter a, is 



related to the symmetry energy C sym , through the rela- 
tion, 



AC sym (Zl Z$\ 

T U Al) [2) 

where, Z\ , A\ and Z 2 , A 2 are the charge and the mass 
numbers from the two systems and T is the temperature. 
This relation provides a simple and straight-forward con- 
nection between the symmetry energy and the fragment 
isotopic yield distribution. 

It must be mentioned that although the above equation 
derived from the statistical and the dynamical models of 
multifragmentation appears similar in form, the physical 
meaning of the terms involved in this equation differ for 
each model. 

1) In statistical models, the Z/A in Eq. (2) corresponds 
to the charge-to-mass ratio of the initial equilibrated frag- 
menting system. Whereas, in dynamical models, it cor- 
responds to the charge-to-mass ratio of the liquid phase 
at a certain time (~ 300 fm/c) during the dynamical evo- 
lution of the colliding systems. 

2) The interpretation of the symmetry energy C sym , in 
dynamical and statistical models also differs significantly. 
The dynamical models relate the symmetry energy in the 
above equation to that of the fragmenting source. The 
statistical models, on the other hand, relate C sym to that 
of the fragments formed at freeze-out. 

These conceptual differences between the statistical 
and the dynamical models arc due to the radically differ- 
ent approaches taken in the interpretation of the multi- 
fragmentation process. The different interpretation has 
also lead to conflicting results from the use of Eq. 2, due 
to the different sequential decay effects predicted for the 
primary fragments by each model. 

The isoscaling parameter a, in Eq. 2 corresponds to 
the hot primary fragments which undergo sequential de- 
cay into cold secondary fragments. These secondary frag- 
ments are the ones that are eventually detected in exper- 
iments. The experimentally determined isoscaling pa- 
rameter must therefore be corrected for the sequential 
decay effect before comparing it to the theoretical mod- 
els. It has been observed that while statistical model 
calculations show no significant change in the isoscaling 
parameter after sequential decay [4fJ, dynamical models 
give contrasting results; with some showing no significant 
changes 1411 . while others showing a change of as much 
as 50% gj 

In this work, we adopt both theoretical approaches 
with their respective interpretations to study the den- 
sity dependence of the symmetry energy. In particular, 
we use the Antisymmetrized Molecular Dynamics (AMD) 
model [3l|, HH and the Statistical Multifragmentation 
Model (SMM) [3(| for this study. A comparison between 
the two can provide useful insight into the physical mean- 
ing of the above equation in the two models. 
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III. EXPERIMENT 
A. Experimental Setup 

The experiments were carried out at the Cyclotron 
Institute of Texas A&M University (TAMU) using the 
K500 Superconducting Cyclotron and the National Su- 
perconducting Cyclotron Laboratory (NSCL) at Michi- 
gan State University (MSU). Targets of 58 Fe (2.3 
mg/cm 2 ) and 58 Ni (1.75 mg/cm 2 ) were bombarded with 
beams of 40 Ar and 40 Ca at 33 and 45 MeV/nucleon for 
the TAMU measurements 44|, and targets of 58 Fe (~ 5 
mg/cm 2 ) and 58 Ni (~ 5 mg/cm 2 ) were bombarded with 
beams of 40 Ar and 40 Ca at 25 and 53 MeV/nucleon for 
the NSCL measurements [iljj . The various combinations 
of target and projectile nuclei allowed for a range of N/Z 
(neutron-to-proton ratio) (1.04 — 1.23) of the system to 
be studied, while keeping the total mass constant (A = 
98). In a separate experiment at TAMU, beams of 58 Ni 
and 58 Fe at 30, 40, and 47 MeV/nucleon were also bom- 
barded on self-supporting 58 Ni and 58 Fe targets. 

The beams in the TAMU measurements were fully 
stripped by allowing them to pass through a thin alu- 
minum foil before being hit at the center of the target 
inside the TAMU Air neutron ball [46|]. Light charged 
particles (Z < 2) and intermediate mass fragments (Z > 
2) were detected using six discrete telescopes placed in- 
side the scattering chamber of the neutron ball at angles 
of 10°, 44°, 72°, 100°, 128° and 148°. Each telescope con- 
sisted of a gas ionization chamber (IC) followed by a pair 
of silicon detectors (Si-Si) and a Csl scintillator detector, 
providing three distinct detector pairs (IC-Si, Si-Si, and 
Si-Csl) for fragment identification. The ionization cham- 
ber was of axial field design and was operated with CF4 
gas at a pressure of 50 Torr. The gaseous medium was 6 
cm thick with a typical threshold of ~ 0.5 MeV/nucleon 
for intermediate mass fragments. The silicon detectors 
had an active area of 5 cm x 5 cm and were each sub- 
divided into four quadrants. The first and second silicon 
detectors in the stack were 0.14 mm and 1 mm thick, 
respectively. The dynamical energy range of the silicon 
pair was ~ 16 - 50 MeV for 4 He and - 90 - 270 MeV 
for 12 C. The Csl scintillator crystals that followed the 
silicon detector pair were 2.54 cm in thickness and were 
read out by photodiodes. Good elemental (Z) identifi- 
cation was achieved for fragments that punched through 
the IC detector and stopped in the first silicon detector. 
Fragments measured in the Si-Si detector pair also had 
good isotopic separation. Fragments that stopped in Csl 
detectors showed isotopic resolution up to Z = 7. The 
trigger for the data acquisition was generated by requir- 
ing a valid hit in one of the silicon detectors. 

The calibration of the IC-Si detectors were carried out 
using the standard alpha sources and by operating the 
IC at various gas pressures. The Si-Si detectors were cal- 
ibrated by measuring the energy deposited by the alpha 
particles in the thin silicon and the punch-through ener- 
gies of different isotopes in the thick silicon. The Si-Csl 



detectors were calibrated by selecting points along the 
different light charged isotopes and determining the en- 
ergy deposited in the Csl crystal from the energy loss in 
the calibrated Si detector. 

The setup for the NSCL experiment consisted of 13 
silicon detector telescopes placed inside the MSU 4ir Ar- 
ray. Four of which were placed at 14°, each of which 
consisted of a 100 /im thick and a 1 mm thick silicon 
surface-barrier detector followed by a 20 cm thick plastic 
scintillator. Five telescopes were placed at 40°, in front 
of the most forward detectors in the main ball of the 4ir 
Array. They each consisted of a 100 /im surface-barrier 
detector followed by a 5 mm lithium drifted silicon de- 
tector. More details can be found in Ref. [45j]. Good 
isotopic resolution was obtained as in TAMU measure- 
ments. 



B. Event Characterization 

The event characterization of the NSCL data was 
accomplished by detection of nearly all the coincident 
charged particles by the MSU 47T Array. Data were ac- 
quired using two different triggers; the bulk of which 
was obtained with the requirement of a valid event in 
one of the silicon telescopes. Additional data were taken 
with a minimum bias 4-7T Array trigger for normalization 
of the event characterization. The impact parameter of 
the event was determined by the mid-rapidity charge de- 
tected in the An Array as discussed in Ref. [13] ■ The 
effectiveness of the centrality cuts was tested by compar- 
ing the multiplicity of events from a minimum bias trigger 
with the multiplicity distribution when a valid fragment 
was detected at 40° 48]. The minimum bias trigger had 
a peak multiplicity of charged particles of one, whereas 
with the requirement of a fragment at 40°, the peak of 
the multiplicity distribution increased to five. 

The event characterization for the TAMU data was 
accomplished by using the 47r neutron ball that sur- 
rounded the detector assembly. The neutron ball con- 
sisted of eleven scintillator tanks segmented in its me- 
dian plane and surrounding the vacuum chamber. The 
upper and the lower tank were 1.5 m diameter hemi- 
spheres. Nine wedge-shaped detectors were sandwiched 
between the hemispheres. All the wedges subtended 40° 
in the horizontal plane. The neutron ball was filled with 
a pseudocumene-based liquid scintillator mixed with 0.3 
% (b.w.) of Gd salt (Gd 2-ethyl hexanoate). Scintilla- 
tions from thermal neutrons captured by Gd were de- 
tected by twenty 5-inch phototubes : five in each hemi- 
sphere, one on each of the identical 40° wedges and two 
on the forward edges. The efficiency with which the neu- 
trons could be detected is about 83%, as measured with 
a 252 Cf source. 

The detected neutrons were used to differentiate be- 
tween the central and peripheral collisions. To un- 
derstand the effectiveness of neutron multiplicity as 
a centrality trigger, simulations were carried out us- 
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ing a hybrid BUU-GEMINI calculations at various im- 
pact parameters for the 40 Ca + 58 Fe reaction at 33 
MeV/nucleon. The simulated neutron multiplicity distri- 
bution was compared with the experimentally measured 
distribution. The multiplicity of neutron for the impact 
parameter b — collisions was found to be higher than 
the b — 5 collision. By gating on the 10% highest neu- 
tron multiplicity events, one could clearly discriminate 
against the peripheral events. 

To determine the contributions from noncentral im- 
pact parameter collisions, neutrons emitted in coinci- 
dence with fragments at 44° and 152° were calculated at 
b = fm and b = 5 fm. The number of events were ad- 
justed for geometrical cross sectional differences. A ratio 
was made between the number of events with a neutron 
multiplicity of at least six, calculated at b = fm, and the 
number of events with the same neutron multiplicity at b 
= 5 fm. The ratios were observed to be 19.0 and 11.1 at 
44° and 1.3 and 2.2 at 152° for 33 and 45 MeV/nucleon 
respectively. At intermediate angles, high neutron multi- 
plicities were observed to be outside the region in which 
6 = 5 fm contributes significantly. At backward angles 
the collisions at b = 5 fm made a larger contribution to 
the neutron multiplicity. 

In addition to the neutron multiplicity distribution, the 
charge distribution of the fragments was also used to in- 
vestigate the contributions from central and mid-impact 
parameter collisions. The 6 = 5 collisions produced es- 
sentially no fragments with charge greater than three in 
the 44° telescope. 

In an earlier work fjij . some analysis of the fragment 
kinetic energy and charge distributions were presented. 
It was shown that at a laboratory angle of 44° the ki- 
netic energy and the charge distributions are well repro- 
duced by the statistical model calculation. Using a mov- 
ing source analysis of the fragment energy spectra, it was 
also shown that the fragments emitted at backward an- 
gles originate from a target-like source, while those emit- 
ted at 44° originate primarily from a composite source. 
In this work, we will concentrate exclusively on data from 
the laboratory angle of 44° , which corresponds to the cen- 
ter of mass angle « 90°, to study the symmetry energy 
and the isoscaling properties of the fragments produced. 
The choice of this angle enables one to select events which 
are predominantly central and undergo bulk multifrag- 
mentation. The contributions to the intermediate mass 
fragments from the projectile-like and target-like sources 
can therefore be assumed to be minimal. 



IV. EXPERIMENTAL RESULTS 
A. Fragment isotopic yield distribution 

The experimentally measured relative isotopic yield 
distributions for the Lithium (left) and Carbon (right) 
elements, in 58 Ni + 58 Ni (star symbols), 58 Ni + 58 Fe 
(square symbols), 58 Fe + 58 Ni (circle symbols) and 58 Fe 
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FIG. 2: Relative yield distribution of the fragments for the 
Lithium (left) and Carbon (right) isotopes in 58 Ni + 58 Ni 
(stars and solid lines), 58 Fe + 58 Ni (circles and dashed lines), 
58 Ni + 58 Fe ( squares an{ j dashed lines), and 58 Fe + 58 Fe (tri- 
angles and dotted lines) reactions at various beam energies. 



+ 58 Fe (triangle symbols) reactions, are shown in Fig. 2 
for beam energies of 30, 40 and 47 MeV/nucleon. Sim- 
ilarly, the isotopic yield distributions for Lithium (left), 
Berillium (center) and Carbon (right) elements, in 40 Ca 
+ 58 Ni (star symbols), 40 Ar + 58 Ni (circle symbols) and 
40 Ar + 58 Fe (square symbols) reactions, are shown in fig- 
ure 3 for beam energies of 25, 33 and 45 MeV/nucleon. 
The isotope distribution for each element in Fig. 3 
shows higher fragment yield for the neutron rich iso- 
topes in 40 Ar + 58 Fe reaction (squares), which has the 
largest neutron-to- proton ratio (N/Z), in comparison to 
the 40 Ca + 58 Ni reaction (stars), which has the small- 
est neutron-to-proton ratio. The yields for the reaction, 
40 Ar + 58 Fe (circles), which has an intermediate value of 
the neutron-to-proton ratio, are in between those of the 
other two reactions. A similar feature is also observed for 
the 58 Ni + 58 Fe, 58 Fe + 58 Ni and 58 Fe + 58 Fe reactions 
shown in Fig. 2. The fragment yield distributions there- 
fore show the isospin dependence of the composite system 
on the fragments produced in the multifragmentation re- 
action. One also observes that the relative difference in 
the yield distribution between the three reactions in each 
figure decreases with increasing beam energy. This is due 
to the secondary de-excitation of the primary fragments, 
a process that becomes important for systems with in- 
creasing neutron-to-proton ratio and excitation energy. 
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FIG. 3: Relative yield distribution of the fragments for 
Lithium (left), Berillium (center) and Carbon (right) isotopes 
in 40 Ca + 58 Ni (stars and solid lines), 40 Ar + 58 Ni (circles and 
dashed lines), and 40 Ar + 58 Fe (squares and dotted lines) re- 
actions at various beam energies. 



B. Isotopic and Isotonic scaling 

As discussed in section II, the ratio of isotope 
yields in two different systems, 1 and 2, i?2i(-/V, Z) = 
Y 2 (N,Z)/Yi(N,Z), follows an exponential dependence 
on the neutron number (A^) and the proton number (Z) 
of the isotopes in relation known as isoscaling. 

In Fig. 4, we show the isotopic yield ratio as a func- 
tion of neutron number N, for Ar + Fe, Ar + Ni and 
Ca + Ni systems at beam energies of 25, 33, 45 and 
53 MeV/nucleon. The left column shows the ratio for 
the 40 Ar + 58 Fc and 40 Ca + 58 Ni pair of reaction and 
the right column shows the ratio for the 40 Ar + 58 Ni and 
40 Ca + 58 Ni pair of reaction. One observes that the ratio 
for each element shows linear behavior in the logarithmic 
plot and aligns with the neighboring element quite well. 
This feature is observed for all the beam energies and 
both pairs of reactions studied. One also observes that 
the alignment of the data points varies with beam ener- 
gies as well as the pairs of reaction. To have a quantita- 
tive estimate of this variation, the ratio for each element 
(Z) was simultaneously fit using an exponential relation 
(shown by the solid lines) to obtain the slope parame- 
ter a. The values of the parameters are shown at the 
top of each panel in the hgurc. The value of the slope 



10 



Ar + Fe / Ca + Ni Ar + Ni / Ca + Ni 



CrT 10 



o 

n: 

> 
l < 

o 

I 10 



1 r 



: a =0.434 

L| 0b — AD W 


a =0.27 : 

cv / IIULIcUll 


■ 1 1 1 1 1 1 1 1 1 1 1 1 1 - 
: a =0.365 

\ ^ 

E* = 33 U 


■ 1 1 1 1 1 1 1 1 1 1 1 1 1 ■ 

a =0.198 ; 

r t^^3^ : 

eV/nucleon 


. i i i i i i i i i i i i i - 
; a = 0.255 : 

F - 4^ V 
Hob " J 


. i i i i i i i i i i i i i . 

a = 0.098 ; 

p\/ /ni iflpon 

C V / 1 1 UL IdUl 1 


■ 1 1 1 1 1 1 1 1 1 1 1 1 1 - 

: a =0.192 

E* = 53 u 


. i i i i i i i i i i i i i : 

a = 0.093 : 

eV/nucleon 

i , , , i , , , i , , , i , , 



2 4 6 



2 4 6 



FIG. 4: Experimental isotopic yield ratios of the fragments as 
a function of neutron number N, for various beam energies. 
The left column correspond to 40 Ar + 58 Fe and 40 Ca + 



Ni 



pair of reactions. The right column correspond to Ar + Ni 
and 40 Ca + 58 Ni pair of reactions. The different symbols 
correspond to Z = 3 (circles), Z — 4 (open stars), Z = 5 
(triangles) , Z — 6 (squares) and Z = 7 (filled stars) elements. 
The lines are the exponential fits to the data as explained in 
the text. 



parameter a is larger for the 40 Ar + 58 Fe and 40 Ca + 
58 Ni reactions, which has a larger difference in the N/Z 
of the systems in the pair, compared to the 40 Ar + 58 Ni 
and 40 Ca + 58 Ni reactions, which has a smaller difference 
in the corresponding N/Z. The a value furthermore de- 
creases with increasing beam energy. A similar feature 
is also observed in Fe + Fe, Fe + Ni, Ni + Fe and Ni 
+ Ni systems. Fig. 5 shows the isotope yield ratios and 
the isotone yield ratios for the Fe + Fe and Ni + Ni re- 
actions for the 30 MeV/nucleon beam energy. A relative 
comparison of how the isoscaling parameter a, evolves as 
a function of beam energy and the isospin of the system 
is shown in Fig. 6. The hgure clearly shows that the a 
value decreases with beam energy from 25 MeV/nucleon 
to 53 MeV/nucleon. In addition, there is also a clear 
drop in the a values with the decrease of the N/Z values 
of the system. 
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FIG. 5: Experimental isotope yield ratios (top) and isotone 
yield ratios (bottom) from 58 Fe + 58 Fe and 58 Ni + 58 Ni re- 
actions as a function of N and Z for 30 MeV/nucleon beam 
energy. The solid lines are fit to the data as discussed in the 
text. 
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Experimental isoscaling parameter a, as a function 
of the beam energy. The solid circles are for the 40 Ar + 58 Fe 
and 40 Ca + 58 Ni reactions. The open triangles are for 58 Fe + 
58 Fe and 58 Ni + 58 Ni reactions. The solid stars are for 40 Ar 
+ 58 Ni and 40 Ca + 58 Ni reactions. The open squares are for 
58 Fc + 58 Ni and 58 Ni + 58 Ni reactions. 



V. THEORETICAL MODEL COMPARISON 

A. Dynamical AMD model 

The Antisymmetrized Molecular Dynamics (AMD) 
[3lL [43l | is a microscopic model that simulates the time 
evolution of a nuclear collision. The colliding system in 
this model is represented in terms of a fully antisym- 
metrized product of Gaussian wave packets. During the 
evolution, the wave packet centroids move according to 
the deterministic equation of motion. The followed state 
of the simulation branches stochastically and successively 
into a huge number of reaction channels. The interactions 
are parameterized in terms of an effective force acting be- 
tween nucleons and the nucleon-nucleon collision cross- 
sections. The advantage of using a dynamical model to 
study the nuclear equation of state is that it allows one 
to understand the functional form of the density depen- 
dence of the symmetry energy at a very fundamental level 
(i.e., from the basic nucleon-nucleon interaction). 

Recently 3l|, the fragment yields from heavy ion col- 
lisions simulated within the Antisymmetrized Molecular 
Dynamics (AMD) calculation were reported to follow a 
scaling behavior of the type shown in Eq. 1 . A linear re- 
lation between the isoscaling parameter a and the differ- 
ence in the isospin asymmetry (Z/A) 2 of the fragments 
(as given in Eq. 2), with appreciably different slopes, 
was predicted for two different forms of the density de- 
pendence of the symmetry energy; a " stiff " dependence 
(obtained from Gogny-AS interaction) and a " soft " de- 
pendence (obtained from Gogny interaction). 

In this section, we compare the experimentally de- 
termined isoscaling parameter with the predictions of 
the AMD model calculation. The isospin asymmetry 
of the fragments for the present systems was estimated 
at t = 300 fm/c of the dynamical evolution using the 
AMD calculation. The values for the fragment asym- 
metry (Z/A) 2 , were obtained by interpolating between 
those calculated for the 40 Ca + 40 Ca, 48 Ca + 48 Ca and 
60 Ca + 60 Ca systems by Ono et al. [3l|. These systems 
are symmetric and nearly similar in charge and mass as 
studied in the present work. Fig. 7 shows the AMD cal- 
culation of the fragment asymmetry, (Z/A) 2 at t = 300 
fm/c, as a function of initial asymmetry at time t = 
fm/c, for two different choices of the nucleon-nucleon in- 
teraction, Gogny and Gogny-AS. The asymmetry values 
for the 40 Ca + 40 Ca, 48 Ca + 48 Ca and 60 Ca + 60 Ca sys- 



tems of Ref. [3l| are shown by solid and hollow square 
symbols for the Gogny and Gogny-AS interaction, re- 
spectively. The lines are the linear fits to the calcula- 
tions. The interpolated values for the present systems 
are shown by the solid circles and triangles, and hollow 
circles and triangles for the Gogny and Gogny-AS inter- 
action, respectively. 

We note that the AMD calculations carried out in Ref. 
[3ll | and shown in Fig. 7 correspond to the beam energy 
of 35 MeV/nucleon. The interpolated values of the asym- 
metries for the present systems obtained from Fig. 7 are 
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FIG. 7: AMD calculations of the fragment asymmetry 
(Z/A) 2 , at t — 300 fm/c for the Gogny (solid line and solid 
squares) and Gogny-AS (dotted line and hollow squares) in- 
teractions at 35 MeV/nucleon. The calculations are taken 
from Ref. [3lJ for the systems shown by the square symbols. 
The lines are linear fit to the square symbols. The other sym- 
bols are the interpolated values for the systems studied in this 
work. 



FIG. 8: Isoscaling parameter a, as a function of the differ- 
ence in fragment asymmetry for 35 MeV/nucleon. The solid 
and the dotted lines are the AMD calculations for the Gogny 
and Gogny-AS interactions, respectively [3l[. The solid and 
the hollow squares, stars, triangles and circles are from the 
present work as described in the text. The other symbols cor- 
responds to data taken from [4!| (asterisks) and [3q ] (crosses, 
diamonds, inverted triangles). 



therefore for the beam energy of 35 MeV/nucleon. In or- 
der to compare the experimentally determined isoscaling 
parameter to that of the calculations, we therefore make 
use of the experimental isoscaling parameter for the beam 
energy of 35 MeV/nucleon using Fig. 6. 

Fig. 8 shows a comparison between the experimentally 
observed a and those from the AMD model calculations 
plotted as a function of the difference in the fragment 
asymmetry for the beam energy of 35 MeV/nucleon. The 
solid and the dotted lines are the AMD model predictions 
for the " soft " (Gogny) and the " stiff " (Gogny-AS) 
form of the density dependence of the symmetry energy, 
respectively. The solid and the hollow symbols (squares, 
stars, triangles and circles) are the results of the present 
study for the two different values of the fragment asym- 
metry, assuming Gogny and Gogny-AS interactions, re- 
spectively. Also shown in the figure are the scaling pa- 
rameters (asterisks, crosses, diamond and inverted tri- 
angle) taken from various other works in the literature 
[39 . l49j . It is observed that the experimentally deter- 
mined a parameter increases linearly with increasing dif- 
ference in the fragment asymmetry of the two systems 
as predicted by the AMD calculation. Also, the data 
points are in closer agreement with those predicted by 
the Gogny-AS interaction (dotted line) than those from 
the usual Gogny force (solid line). 

In the above comparison between the data and the cal- 



culation, the corrections for the isoscaling parameter a 
due to the sequential de-excitation of the fragments are 
not taken into account. The slightly lower values of the 
isoscaling parameters (symbols) from the present mea- 
surements with respect to the Gogny-AS values (dotted 
line) could be due to the small secondary de-excitation 
effect of the fragments not accounted for in this compar- 
ison. Recently, it has been reported by Ono et al. [42| . 
that the sequential decay effect in the dynamical calcu- 
lations can affect the a value by as much as 50 %, and 
the ability to distinguish between the " stiff " and the " 
soft " form of the density dependence of the symmetry 
energy diminishes significantly. The calculations by Ono 
et al, were carried out for the above studied systems 
using the AMD model. However, dynamical calculation 
carried out by Tian et al. [4 ll ]. using Isospin Quantum 
Molecular Dynamic (IQMD) model shows no significant 
difference between the primary and the secondary a. The 
sequential decay effect from the IQMD calculation was 
also carried out for the same systems and beam energy 
as studied by Ono et al. [13] using the AMD model. The 
contrasting results between the two dynamical calcula- 
tions for the same systems and energy currently present 
significant amount of uncertainty in reliably estimating 
the effect of sequential decay using dynamical models. 
One reason for this could be the large discrepancy that 
exists in the determination of the primary fragment exci- 
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tation energy from the current dynamical models. It has 
been shown using another dynamical model (stochastic 
mean field calculation) (see Liu et al. [33[), that it re- 
quires a significantly lower value of the primary fragment 
excitation energy (by as much as 50%), to be able to 
reproduce the experimentally observed fragment isotope 
distribution. 

In the above comparison between the data and the cal- 
culation, we have therefore assumed the effect of the se- 
quential decay to be negligible. A correction of about 10 
- 15 %, as determined and well established from various 
statistical model studies (28[, results in a slight increase 
in the a values bringing them even closer to the dotted 
line. The observed agreement of the experimental data 
with the Gogny-AS type of interaction therefore appears 
to suggest a " stiff " form of the density dependence of 
the symmetry energy. 

Figure 9 shows various forms of the density dependence 
of the symmetry energy in isospin asymmetric nuclear 
matter used by Chen et al. [50| . and those used in the 
present dynamical model analysis. The dot-dashed, dot- 
ted and the dashed curve corresponds to the momentum 
dependent Gogny interactions used by Chen et al., to ex- 
plain the NSCL-MSU isospin diffusion data. Assuming 
that the density dependence of the symmetry energy can 
be parametrized as, 
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FIG. 9: Different forms of the density dependence of the nu- 
clear symmetry energy used in the dynamical analysis of the 
present measurements on isoscaling data and the isospin dif- 
fusion measurements of NSCL-MSU [50(. The curves are as 
described in the text. 
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where C° ym , is the value of the symmetry energy at sat- 
uration density and 7 is the parameter that characterizes 
the stiffness of the symmetry energy, the above depen- 
dences used by Chen et al. can be written as, E s ^ m ~ 
31.6 (p/po) 1 , where, 7 = 1.6, 1.05 and 0.69, respectively. 
The solid curve and the solid point in Fig. 9 correspond 
to those from the Gogny and Gogny-AS interactions used 
to study the isoscaling data in the present work. 

By parameterizing the density dependence of the sym- 
metry energy that explains the present isoscaling data, 
one obtains, C sym (p) « 31.6 (p/p ) 7 , where 7 = 0.69, 
from the dynamical model analysis. 



B. Statistical Multifragmentation Model 

The Statistical Multifragmentation Model (SMM) [HI, 
|52| is based on the assumption of statistical equilibrium 
at a low density freeze-out stage. All breakup channels 
composed of nucleons and excited fragments are taken 
into account and considered as partitions. During each 
partition the conservation of mass, charge, energy, mo- 
mentum and angular momentum is taken into account, 
and the partitions are sampled uniformly in the phase 
space according to their statistical weights using Monte 
Carlo sampling. The Coulomb interaction between the 
fragments is treated in the Wigner-Seitz approximation. 
Light fragments with mass number A < 4 are considered 



as elementary particles with only translational degrees 
of freedom ("nuclear gas"). Fragments with A > 4 are 
treated as heated nuclear liquid drops, and their individ- 
ual free energies Fa,z are parametrized as a sum of the 
volume, surface, Coulomb and symmetry energy. 

For the present study we make use of the SMM ver- 
sion adopted by Botvina et al. [36| . In this version, 
the secondary de-excitation of large fragments with A 
> 16 is described by Weisskopf-type evaporation and 
Bohr-Wheeler- type fission models [SJ, [HI]. The decay 
of smaller fragments is treated with the Fermi-breakup 
model. All ground and nucleon-stable excited states of 
light fragments are taken into account and the popula- 
tion probabilities of these states are calculated according 
to the available phase space (53|. The sequential decay 
effect on the isoscaling parameter in this version of SMM 
has been established to be small and in good agreement 
with other versions of the statistical models. 

Unlike dynamical calculations, the form of the density 
dependence of the symmetry energy is not known a priori, 
but has to be deduced from the systematic correlations 
between the isoscaling parameter, temperature, symme- 
try energy and the density of the multifragmenting sys- 
tem. To build this correlation, we make use of the frag- 
ment yield distributions measured in 58 Ni, 58 Fe + 58 Ni, 
58 Fe reactions at 30, 40 and 47 MeV/nucleon to study 
the isoscaling parameter a, as a function of the excita- 
tion energy of the fragmenting source. The parameter a 
was obtained from the ratio's of the isotopic yields for two 
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different pairs of reactions, 58 Fe + 58 Ni and 58 Ni + 58 Ni, 
and 58 Fe + 58 Fe and 58 Ni + 58 Ni as discussed in section 
4.2. The excitation energy of the source for each beam 
energy was determined by simulating the initial stage of 
the collision dynamics using the Boltzmann-Nordhcim- 
Vlasov (BNV) model calculation The results were 

obtained at a time around 40 - 50 fm/c after the projec- 
tile had fused with the target nuclei and the quadrupole 
moment of the nucleon coordinates (used for identifica- 
tion of the deformation of the system) approached zero. 
These excitation energies were also compared with those 
obtained from the systematic calorimetric measurements 
(see Ref. [55j) for systems with mass (A ~ 100), and sim- 
ilar to those studied in the present work, and are in good 
agreement. Fig. 10 (a) shows the experimental isoscaling 
parameter a, as a function of the excitation energy for 
Fe + Fe and Ni + Ni, and Fe + Ni and Ni + Ni pairs 
of reactions. A systematic decrease in the absolute val- 
ues of the isoscaling parameter with increasing excitation 
energy is observed for both pairs. The a parameters for 
the 58 Fe + 58 Fe and 58 Ni + 58 Ni are about twice as large 
compared to those for the 58 Fe + 58 Ni and 58 Ni + 58 Ni 
pair of reactions. 

The experimental isoscaling parameter was compared 
with the predictions of the Statistical Multifragmentation 
Model (SMM) [UHl] calculations to study their depen- 
dence on the excitation energy and the isospin content. 
The initial parameters such as, the mass, charge and ex- 
citation energy of the fragmenting source for the calcula- 
tion was obtained from the BNV calculations as discussed 
above. The possible uncertainties in the source param- 
eters due to the loss of nucleons during pre-equilibrium 
emission was accounted for by carrying out the calcula- 
tions for smaller source sizes. The break-up density in 
the calculation was taken to be multiplicity-dependent 
and was varied from approximately 1/2 to 1/3 the sat- 
uration density. This was achieved by varying the free 
volume with the excitation energy as shown in Ref. [5l| . 
The form of the dependence was adopted from the work 
of Bondorf et al. [53, HU , (and shown by the solid curve in 
Fig. 10 (d)). It is known that the multiplicity-dependent 
break-up density, which corresponds to a fixed interfrag- 
ment spacing and constant pressure at brea k-up, le ads 
to a pronounced plateau in the caloric curve [53, [58| • A 
constant break-up density would lead to a steeper tem- 
perature versus excitation energy dependence. 

The symmetry energy in the calculation was varied un- 
til a reasonable agreement between the calculated and the 
measured a was obtained. Fig. 10 (a) shows the com- 
parison between the SMM calculated and the measured 
a for both pairs of systems. The dashed curves corre- 
spond to the calculation for the primary fragments and 
the solid curves to the secondary fragments. The width 
in the curve is the measure of the uncertainty in the in- 
puts to the SMM calculation. It is observed that, within 
the given uncertainties, the decrease in the a values with 
increasing excitation energy and decreasing isospin dif- 
ference A(Z/A) 2 , of the systems is well reproduced by 
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FIG. 10: (Color online) Isoscaling parameter a, temperature, 
symmetry energy and density as a function of excitation en- 
ergy for the Fe + Fe and Ni + Ni (inverted triangles), and 
Fe + Ni and Ni + Ni (solid circles) reactions at 30, 40 and 
47 MeV/nucleon. a) Experimental isoscaling parameter as 
a function of excitation energy. The solid and the dashed 
curves are the SMM calculations as discussed in the text, b) 
Temperature as a function of excitation energy. The solid 
stars correspond to the measured values and are taken from 
Ref. [55l ]. The solid and the dashed curve corresponds to the 
Fermi-gas relation. The dotted curve corresponds to the one 
obtained from Eq. 4. c) Symmetry energy as a function of 
excitation energy, d) Density as a function of excitation en- 
ergy. The solid stars correspond to those from Ref. The 
open triangles are those from Ref. [70j |. The solid curve is 
from Ref. 15711. 
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the SMM calculation. One also notes that the effect of 
sequential decay on the isoscaling parameter is small as 
observed in several other studies [40, using statistical 
models. 

We show in Fig. 10 (b), the temperature as a func- 
tion of excitation energy (caloric curve) obtained from 
the above SMM calculation that uses the excitation en- 
ergy dependence of the break-up density to explain the 
observed isoscaling parameters. These are shown by the 
solid and inverted triangle symbols. Also shown in the 
figure are the experimentally measured caloric curve data 
compiled by Natowitz et al. [55} , from various measure- 
ments for this mass range. The data from these mea- 
surements are shown collectively by solid star symbols 
and no distinction is made among them. The Fermi-gas 
model predictions with inverse level density parameter 
K a = 10 (solid and dashed curve), is also shown. It is 
evident from the figure that the temperatures obtained 
from the SMM calculations are in good agreement with 
the overall trend of the caloric curve. Somewhat lower 
value for the temperature is observed when the break-up 
density of the system is kept constant at 1/3 the nor- 
mal nuclear density. By allowing the break-up density 
to evolve with the excitation energy, a near plateau that 
agrees with the experimentally measured caloric curves 
is obtained. This assures that the input parameters used 
in the SMM calculation for comparing with the data are 
reasonable. 

The symmetry energies obtained from the statistical 
model comparison of the experimental isoscaling param- 
eter a, are as shown in Fig. 10 (c). A steady decrease in 
the symmetry energy with increasing excitation energy is 
observed for both pairs of systems. Such a decrease has 
also been observed in several other studies [rjOl. l6lH62l. 1631] . 
We have also estimated the effect of the symmetry en- 
ergy evolving durin g th e sequential de-excitation of the 
primary fragments [60l . [HJ|. These are reflected in the 
large error bars shown in Fig. 10 (c). 

The phase diagram of the multifragmenting system is 
two dimensional and hence the excitation energy depen- 
dence of the temperature (the caloric curve) must take 
into account the density dependence too. Often this de- 
pendence is neglected while studying the caloric curve. 
In the following, we attempt to extract the density of the 
fragmenting system as a function of excitation energy. It 
has been shown by Sobotka et al. [65[, that the plateau 
in the caloric curve could be a consequence of the ther- 
mal expansion of the system at higher excitation energy 
and decreasing density. By assuming that the decrease 
in the breakup density, as taken in the present statisti- 
cal multifragmentation calculation, can be approximated 
by the expanding Fermi gas model, and furthermore the 
temperature in Eq. 2 and the temperature in the Fermi- 
gas relation are related, one can extract the density as a 
function of excitation energy using the relation 

T = yjK (p/ Po )V3E* (4) 



> 50 




FIG. 11: (Color online) Symmetry energy as a function of 
density for the Fe + Fe and Ni + Ni pair of reaction (inverted 
triangles), and Fe + Ni and Ni + Ni pair of reactions (solid 
circles) for the 30, 40 and 47 MeV/nucleon. The solid curve is 
the dependence obtained form the dynamical model analysis 
as explained in the text. 



In the above expression, the momentum and the fre- 
quency dependent factors in the effective mass ratio are 
taken to be one as expected at high excitation energies 
and low densities studied in this work [6a , 1671 . [68| . 

The resulting densities for the two pairs of systems are 
shown in Fig. 10 (d) by the solid circles and inverted 
triangles. For comparison, we also show the break-up 
densities obtained from the analysis of the apparent level 
density parameters required to fit the measured caloric 
curve by Natowitz et al. (6|| , and those obtained by Viola 
et al. 70j from the Coulomb barrier systematics that are 
required to fit the measured intermediate mass fragment 
kinetic energy spectra. One observes that the present re- 
sults obtained by requiring to fit the measured isoscaling 
parameters and the caloric curve are in good agreement 
with those obtained by Natowitz et al. The figure also 
shows the fixed freeze-out density of 1/3 (dashed line) 
and 1/6 (dotted line) of the saturation density assumed in 
various statistical model comparisons. The caloric curve 
obtained using the above densities and excitation ener- 
gies (shown by solid stars, circles and the triangles) with 
K Q = 10 in Eq. 4, is shown by the dotted curve in Fig. 
10 (b). The small discrepancy between the dotted curve 
and the data (solid stars) below 4 MeV/nucleon is due 
to the approximate nature of Eq. 4 being used. 

It is therefore evident from figure 10 (a), (b), (c) and 
(d) that the decrease in the experimental isoscaling pa- 
rameter a, symmetry energy, break-up density, and the 
flattening of the temperature with increasing excitation 
energy are all correlated. One can thus conclude that the 
expansion of the system during the multifragmentation 
process leads to a decrease in the isoscaling parameter, 
decrease in the symmetry energy and density, and the 
flattening of the caloric curve. 
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TABLE I: Parameterized form of the density dependence of the symmetry energy obtained from various independent studies. 



Reference 


Parametrization 


Studies 


Fuchs et al. [85] 


32.9(p/p D ) u 9 


Relativistic Dirac-Brueckner calculation 


Heiselberg et al. [82] 


32.0(p/p o ) a6 ° 


Variational calculation 


Danielewicz et al. [81] 


31(33)(p/p o ) ' j5lu - 7J) 


BE, skin, isospin analog states 


Tsang et al. [79] 


12.125(p/p D ) 2 


Isospin diffusion 


ntipn et nl r^fll 

V_JllCll of/ lit. 


Ol.O\p/ Po) 


Tc/~»qt\iti nifriicnrin 


Li et al. [80] 


31.6(p/p o ) 0M 


Isospin diffusion 


Piekarewicz aZ. [77, 78] 


32.7(p/p ) ' 64 


Giant resonances 


Shetty et aZ. [73, 74, 75] 


31.6(p/p ) ' 69 


Isotopic distribution 


Famiano et al. [87] 


32.0(p/p o ) 055 


neutron-proton emission ratio 


Tsang et aZ. [28] 


23.4(p/p ) 06 


Isotopic distribution 



From the above correlation between the symmetry en- 
ergy as a function of excitation energy, and the density 
as a function of excitation energy, we obtain the sym- 
metry energy as a function of density. This is shown by 
the inverted triangles and solid circles in Fig. 11 for the 
Fe + Fe and Ni + Ni, and the Fe + Ni and Ni + Ni 
pair of reactions. The temperature in the present work 
remains nearly constant for the range of excitation ener- 
gies studied, the observed decrease in the symmetry en- 
ergy with increasing excitation energy is therefore a con- 
sequence of decreasing density. This is also supported 
by microscopic calculations which shows an extremely 
slow evolution of the symmetry energy with temperature 
[Til . [72[. The evolution is practically negligible for the 
temperature range studied in this work. The solid curve 
in Fig. 11 corresponds to the dependence C sym (p) — 
31.6 (p/p ) ' 69 MeV, obtained from the dynamical Anti- 
symmetrized Molecular Dynamic (AMD) calculation, as 
discussed in the previous section. It is thus observed that 
the dynamical and statistical models lead to similar den- 
sity dependence of the symmetry energy. 



VI. COMPARISON WITH OTHER 
INDEPENDENT STUDIES 

In the following, we compare the form of the den- 
sity dependence of the symmetry energy obtained from 
the present experimentally measured isoscaling param- 
eter using the statistical and the dynamical multifrag- 
mentation models with several other recent independent 
studies. Fig. 12 shows this comparison. The green solid 
curve corresponds to the one obtained from Gogny-AS 
interaction in d yna mical AMD model that explains the 
present results |73l . |74| . assuming the sequential decay 
effect to be small. The inverted triangle and the circle 
symbols also correspond to the present measurements ob- 
tained by c omp aring with the statistical multifragmenta- 
tion model [75j | . The red dashed curve corresponds to the 
one obtained recently from an accurately calibrated rela- 
tivistic mean field interaction, used for describing the Gi- 



ant Monopole Resonance (GMR) in 90 Zr and 208 Pb, and 
the IVGDR in 208 Pb by Piekarewicz et al. [zl Izl- 
The pink dot-dashed curve correspond to the one used 
to explain the isospin diffusion results of NSCL-MSU us- 
ing the isospin dependent Boltzmann-Uehling-Uhlenbeck 
(IBUU) model by Tsang et al. [eJ. The blue dot-dashed 
curve also corresponds to the one used for explaining the 
isospin diffusion data of NSCL-MSU by Chen et al. [U, 
but with the momentum dependence of the interaction 
included in the IBUU calculation. This dependence has 
been further modified to include the isospin dependence 
of the in-medium nucleon-nucleon cross-section by Li et 
al. [80( , and is in good agreement with the present study. 
The shaded region in the figure corresponds to those ob- 
tained by constraining the binding energy, neutron skin 
thickness and isospin analogue state in finite nuclei us- 
ing the mass formula of Danielewicz [8l|. The yellow 
solid curve correspond to the parameterization adopted 
by Heiselberg et al. [12] in their studies on neutron stars. 
By fitting earlier predictions of the variational calcula- 
tions by Akmal et al. 83, 84J, where the many-body and 
special relativistic corrections are progressively incorpo- 
rated, Heiselberg and Hjorth-Jensen obtained a value of 
C° ym — 32 MeV and 7 = 0.6, similar to those obtained 
from the present measurements. A similar result is also 
obtained from the relativistic Dirac-Brueckner calcula- 
tion, with C° sym = 32.9 MeV and 7 = 0.59 [85]. The 
Dirac-Brueckner is an "ab-initio" calculation based on 
nucleon-nucleon interaction with Bonn A type potential 
instead of the AV18 potential used in the variational cal- 
culation of Ref. [HI . The density dependence of the sym- 
metry energy has also been studied in the framework of 
expanding emitting source (EES) model by Tsang et al. 
[281 ] . where a power law dependence of the form G sym (p) 
= 23.4(p/p ) 7 , with 7 = 0.6 was obtained. This depen- 
dence (shown by the black dotted curve) is significantly 
softer than other dependences shown in the figure. The 
solid square point in the figure correspond to the value 
of symmetry energy obtained by fitting the experimental 
differential cross-section data in a charge exchange reac- 
tion using the isospin dependent CDM3Y6 interaction of 
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FIG. 12: (Color online) Comparison between the results on the density dependence of the symmetry energy obtained from 
various different studies. The various curves and the symbols are described in the text. 



the optical potential by Khoa et al. 86]. 

An alternate observable, the double neutron/proton 
ratio of nucleons taken from two reaction systems using 
four isotopes of the same element, has recently been pro- 
posed as a probe to study the density dependence of the 
symmetry energy [87| • This observable is expected to be 
more robust than the isoscaling observable. It was shown 
recently [87| that the experimentally determined double- 
ratio for the 124 Sn + 124 Sn reaction to that for the 112 Sn 
+ 112 Sn reaction, results in a dependence with 7 = 0.5 
(shown by black dashed curve in Fig. 12), when com- 
pared to the predictions of the IBUU transport model 
calculations. This observation is in close agreement with 
other studies discussed above. However, this dependence 
has been obtained by using the momentum independent 
calculation of Ref. [88]. A more recent calculation [89[ 
using a BUU transport model that includes momentum 
dependent interaction show significantly lower values for 
the double neutron/proton ratio of free nucleons com- 
pared to the one reported by Famiano et al. 

The parameterized forms of the density dependence of 
the symmetry energy obtained from all the above men- 
tioned studies are as given in Table I. The close agree- 
ment between various independent studies show that a 
constraint on the density dependence of the symmetry 
energy, given as C sym (p) = C° sym (p/ ' p )~< , where C° sym re 



31 - 33 MeV and 7 re 0.55 - 0.69 can be obtained. 

VII. DISCUSSION 

We make the following observations from the above 
comparison between the statistical and the dynamical 
model analysis : 

1) Assuming a negligibly small sequential decay effect, 
the form of the density dependence of the symmetry en- 
ergy obtained from the dynamical model analysis is in 
good agreement with the one obtained from the statistical 
model analysis: As mentioned earlier, the sequential de- 
cay effect among various dynamical model calculations is 
still a subject of debate. The statistical models however 
consistently show small sequential decay effect. If the se- 
quential decay in both the dynamical and the statistical 
model is determined by the excitation energy, charge (Z) 
and mass (A) of the fragments, and not by the process 
that leads to these fragments, the de-excitation of the 
fragments must lead to a same amount of change in the 
isoscaling parameter (either a large change or no change 
at all). It is therefore unrealistic to assume that the se- 
quential decay effect is different in the dynamical and the 
statistical model calculations. One comparison by Hudan 
et al. 9C| , show good agreement between the experimen- 
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tally determined primary fragment excitation energy and 
those calculated using the statistical multifragmentation 
model (see table II of Ref. [9(|)- Furthermore, if dy- 
namical and statistical models are merely two different 
ways of interpreting the same multifragmentation process 
(i.e., one simulating the entire process from the forma- 
tion to the breakup stage, and the other simulating only 
the later breakup stage), the isoscaling parameter from 
both interpretation must lead to consistent results. It 
is well known, and as discussed in section II, that both 
interpretations predict isoscaling in multifragmentation. 
As discussed in section V A, the apparent disagreement 
between the sequential decay effect in statistical and dy- 
namical models, could be due to the large discrepancy 
that exists in the determination of the primary fragment 
excitation energy from current dynamical models. 

It has been argued [9l[ that the effect of sequential 
decay on the isoscaling parameter a, in statistical multi- 
fragmentation model depends not only on the excitation 
energy but also on the value of the symmetry energy. The 
fragments in their primary stage are usually hot and the 
properties of hot nuclei (i.e., their binding energy and 
mass) differ from those of cold nuclei. If hot fragments 
in the freeze-out configuration have smaller symmetry 
energy, their mass at the beginning of the sequential de- 
excitation will be different and this effect should be taken 
into account. At smaller values of the symmetry energy 
the sequential decay effect can be large. In order to study 
the effect of symmetry energy evolution on isoscaling pa- 
rameter during sequential decay, we have adopted in this 
work a phenomenological approach of Buyukcizmeci et 
al. [Hj]. Fig. 13 and 14 shows the primary and the sec- 
ondary isoscaling parameter as a function of symmetry 
energy calculated from the statistical multifragmentation 
model (SMM) for the Ar + Ni and Ca + Ni, and Ar + Fe 
and Ca + Ni pair of reactions. The various panels from 
top to bottom correspond to different system excitation 
energies. Fig. 13 shows the result of the calculations 
where the symmetry energy is kept fixed, and Fig. 14 
shows the result for the calculations where the symme- 
try energy is varied during the de-excitation process. The 
dashed lines in the figure correspond to the primary frag- 
ments (Eq. 2) and the solid lines to the secondary frag- 
ments. It is observed that there is no significant change 
in the primary and the secondary alpha. 

2) The result of the statistical model analysis is in good 
agreement with other independent studies : A comparison 
between the density dependence of the symmetry energy 
obtained from the statistical model analysis (for which 
the sequential decay effect is known to be small) and 
other independent studies shows excellent agreement. 

3) The isoscaling parameter probes the property of in- 
finite nuclear matter : The symmetry energy obtained 
from dynamical model analysis (shown by the solid curve 
in Fig. 11) relates to the volume part of the symmetry en- 
ergy as in infinite nuclear matter, whereas, the symmetry 
energy obtained from the statistical model analysis (solid 
circles and inverted triangles in Fig. 11) relates to the 
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FIG. 13: SMM calculated isoscaling parameter a as a func- 
tion of symmetry energy for various excitation energies. The 
open circles joined by dotted lines correspond to the primary 
fragments and the open stars joined by solid lines to the sec- 
ondary fragments. The left column shows the calculation for 
40 Ar + 58 Ni and 40 Ca + 58 Ni pair, and the right column for 
the 40 Ar + 58 Fe and 40 Ca + 58 Ni pair. 



fragment that is finite and has surface contribution. The 
similarity between the two can probably be understood 
in terms of the weakening of the surface symmetry free 
energy when the fragments are being formed. During 
the density fluctuation in uniform low density matter, 
the fragments are not completely isolated and continue 
to interact with each other, resulting in a decrease in 
the surface contribution as predicted by Ono et al. [43| . 
It must be mentioned that by fitting the nuclear masses 
with mass formula, a volume contribution to the symme- 
try energy of about 31-33 MeV and surface contribution 
of about 11-13 MeV was obtained by Danielewicz [92j 
for nuclei at normal density. Using the constraint ob- 
tained for the volume part of the symmetry energy from 
the present study, and following the expression for the 
symmetry energy of finite nuclei by Danielewicz, we write 
the general expression for the density dependence of the 
symmetry energy as, 



Sa(p) = 



a(p/p ) 7 



1 + Hp/po 



(5) 



where, a = C° sym = 31-33 MeV, 7 = 0.55 - 0.69 and 
a//3 = 2.6 - 3.0. The quantities a and (3 are the volume 
and the surface symmetry energy at normal nuclear den- 
sity. The above equation reduces to Eq. 3 for infinite 
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nuclear matter in the limit of A —> oo, and to the sym- 
metry energy of finite nuclei for p = p a . The ratio of the 
volume symmetry energy to the surface symmetry energy 
(ct/f3) is closely related to the neutron skin thickness. De- 
pending upon how the nuclear surface and the Coulomb 
contribution is treated, two different correlations between 
the volume and the surface symmetry energy have been 
predicted from fits to nuclear masses. Experimental 
masses and neutron skin thickness measurements for nu- 
clei with N/Z > 1 should provide tighter constraint on 
the surface-volume correlation. 

4) The density dependence of the symmetry energy ob- 
tained using the statistical model approach is consistent 
with other experimentally determined observables : In the 
past, attempts have been made to study the density de- 
pendence of the symmetry energy by looking at specific 
observables and comparing them with the predictions of 
the dynamical models. Such an approach attempts to ex- 
plain the observable of interest without trying to simul- 
taneously explain other properties, such as, the temper- 
ature, density and excitation of the fragmenting system. 
This has lead to a variety of different dependences with- 
out any clue to what density is being probed. It might not 
be straightforward to distinguish different realistic EOS 
interactions using dynamical models, due to the large 
uncertainties that currently exist in the sequential decay 
effects for these models. But the idea of extracting infor- 
mation on the symmetry energy from the point of view 
of the basic nucleon-nucleon interaction is a very pow- 



erful approach. On the other hand, the determination 
of the density dependence of the symmetry energy from 
statistical model analysis by simultaneously explaining 
the isoscaling parameter, caloric curve and the density 
as a function of excitation energy is a reverse approach. 
This approach attempts to explain the experimental ob- 
servables without any prior knowledge of the governing 
interaction and arrives at a dependence which can then 
be compared with those predicted from the basic inter- 
actions. 

5) Symmetry energy determined from the multifrag- 
mentation study is lower than that of normal nuclei 
: Theoretical many-body calculations [3, [H, [94], [9f| 
and those from the empirical liquid-drop mass formula 
[96l l97j predict symmetry energy near normal nuclear 
density and temperature to be around 30 - 32 MeV. As- 
suming a negligible evolution of the symmetry energy as 
a function of temperature, as shown in Ref. [Zl[Z3, the 
present statistical model analysis yields symmetry energy 
of the order of 18 - 20 MeV at half the normal nuclear 
density. 

6) The above constraint on the density dependence of 
the symmetry energy has important implications for as- 
trophysical and nuclear physics studies : 

a) Neutron skin thickness: It has been shown [98j that 
an empirical fit to a large number of mean-field calcu- 
lations yield neutron skin thickness for 208 Pb nucleus, 
R n - R p ~ (0.227 + 0.06) fin, where 7 is the exponent 
that determines the stiffness of the density dependence 
of the symmetry energy. From the above comparison, as- 
suming only those density dependences of the symmetry 
energy which have symmetry energy at normal nuclear 
density 31 - 33 MeV, one obtains a neutron skin thick- 
ness of 0.18 - 0.21 fin. An accurate determination of 
the neutron skin thickness from the parity-violating elec- 
tron scattering measurement [18| will provide a unique 
observational constraint on the thickness of the neutron 
skin of a heavy nucleus. The above constraint also leads 
to symmetric matter nuclear compressibility of K ~ 230 
MeV. 

b) Neutron star mass and radius : The constraint also 
predicts a limiting neutron star mass of M max = 1.72 
solar mass and a radius, R = 11 - 13 km for the " canon- 
ical " neutron star. Recent observations of pulsar-white 
dwarf binaries at the Arecibo observatory suggest a pul- 



2.1 



+0.4 
-0.5 



solar mass 



sar mass for PSRJ0751+1807 of M 
at a 95% confidence level [99(. 

c) Neutron star cooling : Furthermore, it predicts a 
direct URCA cooling for neutron stars above 1.4 times 
the solar mass. If such is the case, then the enhanced 
cooling of a M = 1.4 solar mass neutron star may provide 
strong evidence in favor of exotic matter in the core of a 
neutron star. 

These results have important implications for nuclear 
astrophysics and future experiments probing the prop- 
erties of nuclei using beams of neutron-rich nuclei. The 
above constraint was obtained by studying the low den- 
sity behavior of nuclear matter. Measurements at den- 
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sities above normal nuclear matter should further con- 
straint the form of the symmetry energy. Such measure- 
ments should yield consistent results when extrapolated 
to low densities. 



VIII. SUMMARY AND CONCLUSIONS 

In summary, a number of reactions were studied to in- 
vestigate the density dependence of the symmetry energy 
in the equation of state of isospin asymmetric nuclear 
matter. The results were analyzed within the framework 
of the dynamical and the statistical models of multifrag- 
mentation. It is observed that a dependence of the form 
C sym (p) = 31.6 (p/po) ' 69 MeV, agrees reasonably with 
the experimental data indicating that a " stiff " form of 
the symmetry energy provides a better description of the 



nuclear matter EOS at sub-nuclear densities. A compar- 
ison with several other independent studies shows that 
a constraint on the density dependence of the symmetry 
energy given as C sym (p) = C° sym {p/ p )~< , where C° sym w 
31 - 33 MeV and 7 0.55 - 0.69, can thereby be obtained. 
The present observation has important implications for 
astrophysics, as well as, nuclear physics studies to be car- 
ried out at future radioactive beam facilities worldwide. 
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